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Abstract. We develop a General Fluctuation Formula for phase variables that are odd under 
time reversal. Simulations are used to verify the new formula. 



INTRODUCTION 

The Fluctuation Theorem (FT) (1-10) gives a simple relationship for the probability 
of observing trajectory segments with a specific time averaged value of a given phase 
function, to the probability of observing segments with a negative value: that is, for a 
phase function ^{T), it gives an exact expression for Pr((|)(x) = A) / Pr((|)(x) = -A), 
where (|)(x) is the phase function averaged over a time x. Thus far, FT has only been 
obtained for a restricted set of phase functions (7) and originally work focussed on the 
isoenergetic dynamics of a microcanonical ensemble of initial phaseiT. The only phase 
function that was considered was the time-averaged rate of entropy production per unit 
volume, a(x) which is directly proportional to the dissipative flux and the phase space 
contraction rate (1-7). This fluctuation formula is of particular interest because it gives an 
analytic expression of the probability that, for a finite system and for a finite time, the 
Second Law of thermodynamics will be violated. In addition, FT has been derived from 
the Sinai-Ruelle-Bowen measure (5,6) and hence provides support for the applicability 
of this measure to nonequilibrium steady state systems. More recently, FT for other 
ensembles have been obtained (7). In these cases FT was not expressed in terms of the 
phase space contraction but only in terms of the dissipative flux and an analytic 
expression is obtained which is consistent with the Second Law of thermodynamics. In 
the present paper we obtain General Fluctuation Formulae for any phase variable that is 
odd under time-reversal symmetry. 



DERIVATION OF THE GENERALISED FLUCTUATION 

THEOREM (GET) 

Consider an N-particle system with coordinates and peculiar momenta, 
{qi,q2,..qN,Pi,"PN} -(q^P) -T. The internal energy of the system is 

N 

H„ = V pf / 2m -I- 4>(q) = K -I- $ where ^(q) is the interparticle potential energy, which 

i=i 



is a function of the coordinates, q, of all the particles and K is the total peculiar kinetic 
energy. In the presence of an external field Fg, the thermostatted equations of motion are 
taken to be, 

q-p,/m + Q(r)F, 

Pi = F,(q) + Di(r)F^-a(r)p, 

where Fi(q) = -3$(q)/3q; and a is the thermostat multiplier which in this case is 
applied to the peculiar momenta, and Ci and Di represent the coupling of the system to 
the field. The dissipative flux J, is defined in terms of the adiabatic time derivative of the 
internal energy, H^'^ = - J(r)VF^ . 

Let (|)(r) be an arbitrary phase function and define the time average 

(^.(x) = -rds^(r,(s)). (2) 

At t = the phase space volume occupied by a contiguous bundle of trajectories for 
which -{rjA < (|)i(x) < A + 5A[ is given by 5V(r(0),0) and at time x these phase points 

will occupy a volume 5V(r(x),x) = 5V(r(0),0)e'^'^^'^ where A(x) is the time-averaged 
phase space compression factor along these trajectories. We denote (|)(x) =< <^^(t) >,pj , 

that is the average value of <^^(t) over the set of contiguous trajectories, {F;}. 

If the dynamics is reversible, there will be a contiguous set of initial phases |r*(0)|, 
given by r*(0) = M^(r,(x)), that will occupy a volume 5V(r*(0),0) = 5V(r(x),x) = 
5V(r(0),0)e'^*^'^ along which the time-averaged value of the phase function is 
{|)i*(x) = M ((|)i(x)) . For any ^^(T) that is odd under time reversal, (|)i*(x) = -(|)j(x). 

The probability ratio of observing trajectories originating in an initial phase volume 
and its conjugate phase volume will be related to the initial phase space distribution 
function and the size of the volume elements: 

Pr(5v(r(0),o)) f(r(0),o)5v(r(0),o) 



Pr(5v(r (0),o)) f(r (0),o)5v(r (0),o) 



f(r(0),o)^_^(,). 



f(r(x),o) 



where 



_ gB(i:)x 



B(x)x.lnff^^^«)'«)^ 



f(r(x),o). 



j\(r(s))ds, (4) 



then ln(f(r(0),0)/f(r(x),0)) = 0, and B(x)x = -A(r(x))x = -(3J(x)VF^x where 
(3 = dN / 2K and d is the Cartesian dimensionality. If the initial ensemble is canonical, 
then B(x)x = |3($(x) - $(0)) - A(x)x = -|3J(x)VF^x . 

It is possible that there are non-contiguous bundles of trajectories for which 
■{r; A < (|)^j < A + SAj-, and since these bundles may have different values B(x) the 
probability ratio (Eqn (3)) may differ for each bundle. The probability of observing a 
trajectory for which A < (|)(x) < A + 5A, is obtained by summing over the probabilities 
of observing these m=l,M non-contiguous volume elements, 5V^(r(0),0). If the phase 
function is odd under time reversal symmetry, i then the ratio of the probability of 
observing trajectories for which A < (|)(x) < A -I- 5A to the probability of observing 
conjugate trajectories, for which -A < (|)(x) < -A -I- 5A is, 

M 

y Pr(5V (r(0),0) 
Pr((^(x) = A) _ ;^, -"^ ^ ^' ^ 



Pr((^(x) = -A) £pr(5vjr(0),0)) 



m=l 



M 

J^Pr(5Vjr(0),0)) 



m=l 



£pr(5Vjr(0),0)e-^(^'^ 



(5) 



g-5(x)x\-' 

/(t)(T) = A 

where the notation {...)-t (^wa refers to the ensemble average over contiguous trajectory 

bundles for which (|)(x) = A . Equation (5) gives the ratio of the measure of those phase 
space trajectories for which (|)(x) = A to the measure of those trajectories for which 
(|)(x) = -A. This is the Generalised Fluctuation Formula for any phase variable (|)(x) 
that is odd under time reversal. The actual form of (|)(x) is quite arbitrary. 
For isoenergetic dynamics initiated from a microcanonical ensemble (7), 

Pr(^(x) = A) _/^A(x)xr^ ^/ vf,pI(x)t\-i ^g^ 



Pr((^(x) = -A) \ /<^(z)=A \ /(t)(T)=A 

likewise for isokinetic or Nose-Hoover dynamics initiated from a canonical ensemble 

(7), 



1) If the phase variable is even, the we obtain the trivial relationship 

-1 _ Pr((^(T) = A) 
i(x)=A Pr(^(T) = A) 



^_p(,),v-i _ Pr(yT) = A) _ ^ 



-1 



^/gVFjJWx, _ (7^ 



Pr(4)(T) = A) 
Pr((^(x) = -A) \ /$(t)=A 

Formulas for other ergodically consistent ensembles can be obtained in a similar manner 

(7). 



CORRELATIONS BETWEEN TIME AVERAGES AND THE 

ENTROPY PRODUCTION 

Equation (5) is the generalised fluctuation formula for properties, (|)(x) which are 
odd under time reversal mapping. The time averaged entropy production rate per unit 

volume, a(x) is also odd under time reversal symmetry, however (|)(x) does not 

necessarily change sign with o(x) . For example, the kinetic component of the shear 

stress, Pxy(x) is odd under time reversal, but does not necessarily have the same sign as 

o(x). One has to include the possibility that when (|)(x)<0, a(x) > and when (|)(x) 
>0, a(x) < . The probability ratio in Equation (5) can be written as 

Pr((^(x) = A) Pr(a>0)Pr((^(x) = A;o>0) + Pr(a<0)Pr((^(x) = A;o<0) ,^, 

(o) 



Pr(4)(x) = -A) Pr(a > 0)Pr(4)(x) = -A;a > 0) + Pr(o < 0)Pr(4)(x) = -A;o < 0) 

We can consider Equation (8) in more detail. The probability of observing A, 
Pr((|)(x) = A) has been decomposed into two terms. The first "direct" term is the product 
of Pr(o(x)>0), the probability of observing a trajectory with positive entropy 
production, and Pr((|)(x) = A;o(x) > 0), the probability of observing (|)(x) = A on a 
segment with a(x) > . The second term can be considered a "cross term" and accounts 
for cases when o(x) < and (|)(x)>0. Likewise, Pr((|)(x) = -A) is expressed in terms 
of a direct term Pr(a(x) < 0)Pr((|)(x) = -A;a(x) < 0) which accounts for trajectories 
with (^(x)<0, o(x) < 0, and a "cross term" Pr(o(x) > 0)Pr((^(x) = -A;o(x) > 0). We 
can describe the degree of correlation between (|)(x) and o(x) by considering the direct 
and cross term contributions in the probability ratio. Choices of (|)(x) that are highly 
correlated with o(x) imply that the actual value of a(x) strongly depends on A; in the 
limit when (|)(x) = o(x) , the cross terms are exactly zero and 

Pr((^(x) = A) _ Pr(a(x)>0)Pr((^(x) = A;a(x)>0) 



Pr(4)(x) = -A) Pr(o(x)<0)Pr(4)(x) = -A;o(x)<0). (9) 

= exp[-A(x)x] 

In cases where (|)(x) is uncorrelated with o(x), Pr((|)(x) = A;o(x) > 0) = 
Pr((^(x) = A;a(x)<0) and Pr((^(x) = -A;o(x)> 0) = Pr((^(x) = -A;o(x) <0), thus 



Pr(j(x) = A) /^A(.)x\-^ 
Pr(4)(x) = -A) \ /$(T)=A 



where (e '^'''^''■\ is the Kawasaki normalization factor (3). 

(^(x)=A 



RESULTS 

To demonstrate the validity of GFT we have selected five possible phase functions 
(|)(x) that are odd under time reversal. Different components related to the shear stress, a 
non-linear function of the shear stress, and a test function designed to be completely 
unrelated to the entropy production will be considered. Numerical testing of GFT was 
done with NEMD simulations of N=32 WCA particles in two Cartesian directions under 
isoenergetic shear flow using the SLLOD equations of motion for planar Couette flow 
(11). The equations of motion are 

qi=Pi/m + i7yi 
Pi=Fi-i7Pyi-api' 

and Y is the applied strain rate and a is the Gaussian isoenergetic ergostat multiplier. In 
all simulations the number density was n=0.8, the internal energy per particle was 
E/N=1.056 , Y=0.2, and time averaging interval was x=0.08. The isoenergetic SLLOD 
equations of motion are reversible. If we choose the strain rate to be invariant under the 
time reversal mapping, the dissipative flux and the entropy production are odd and the 
reverse trajectory is generated by carrying out the Kawasaki mapping of a point in phase 
space at time t, MK(x, y, px, Py) = (x, -y, -px, Py) and then propagating forward in time 
(2,11). 

As previously discussed, general phase functions that are antisymmetric with respect 
to time reversal need not be correlated with the entropy production. To demonstrate this, 
we consider four functions with varying degrees of entropy production correlation, 

beginning with a highly correlated function ( PP^^^ (x) itself), then strongly correlated (the 

~N 

potential component of [3P (x), Y^PxiC'^)), weakly correlated (the kinetic component of 



PPxy('^)' P^xy^'^)) and finally an uncorrelated^ phase function (^Pxi('^), where Pxi is 



the x-component of the peculiar momenta for the i particle) 
The GFT for H^) = pPxy(T^) can be written as 



2) In practice, it is difficult to choose a "perfectly uncorrelated" phase function, therefore strickly 
speaking, by "uncorrelated" we mean that there is little or no measureable correlation between the 
phase function, and the time averaged entropy production. 



Pr(|3P^y(x) = A) _/ pp^^(x)VF^. 



;K 



Pr(|3P^ rx) = -A) 



xy 



pP,,(t)=A 



J|3P^,(X)+PP'^XT)]VF„X^ 



pP.,(t)=A 



(12) 



AVFx/ (3P,r'c)VF„T 

e Me' 



pPx,(^)=A 



P K 

where the correlation between [3P (x) and (3P (x) determines the "degree" of 
correlation between |3P (x) and o(x). A similar formula can be constructed for 

Pr(pP^y(x) = A)/Pr(pP^y(x) = -A). 

— p — p 

In Figure 1 we show that Pr((3P (x) = A)/Pr((3P (x) = -A) agrees with 



BP^/'OVF^T 



-1 



pP,,(t)=A 



, as predicted by the GFT (Eqn (5)). 
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FIGURE 1. A plot of Pr(PP^y(T) = A) / Pr(PP^y(T) = -A) versus A for a system of 32 WCA 
particles undergoing shear flow. The crosses correspond to the directly determined probability ratio, and 

^|3P„VF„t^"' 



the open circles correspond to ( e 



pP„(t)=A 



The inset shows the relationship between PP (x) 



xy-' 



and PP (T)for a randomly selected set of transient trajectories, giving an indication of the relative direct 



and cross contributions to the probability ratio. The line shows the results expected if pP (T) = PP (x) 
(perfect correlation). 



The observed linear slope is also very close to the slope found for FT using 

Pxy(x) = A (shown as the solid line), indicating that the potential contribution of the 

shear stress is highly correlated to the entropy production rate. This is hardly surprising 

since at this thermodynamic state point 90% of the shear stress comes from the potential 

— — p 

rather than kinetic, contributions. In the inset we show [3P (x) and |3P (x) for a 

randomly selected set of time averaged transient trajectories. Those points within regions 
I and II correspond to "direct" contributions to the probability ratio, that is those 

— p _ — p _ 

trajectories with A >0, |3P (x) = A, a(x)>0, and |3P (x) = -A, a(x)<0 

respectively. We see that in this highly correlated example, almost all the trajectories 

— p _ 

reside in these regions, whereas the "cross" contributions, (3P (x) = A, a(x) < and 



xy' 



(3P (x) = -A, o(x) > are near zero. 



;K 



Figure 2 shows a similar plot for (3P (x). Again GFT (Eqn (5)) is verified although 
the scatter is much larger. If one incorrectly assumed that the constrained average. 



e'f'^'-'^^^'^-^V , could be replaced by (e'^^'^'^^'^^"^^ 



e^''"'' the result 



would be that given by the line. Clearly, this is incorrect. An examination of a similar 

randomly selected set of transient trajectories as shown in the inset indicates that a large 

proportion of the trajectories originate from the "cross" contributions of region II where 

— P 

pPj^ (x) = A , a(x) < , clearly demonstrating that although the kinetic contribution to 
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FIGURE 2. A plot of PrCPP^x) = A) / Pr(PP^ (x) : 



-A) versus A for a system of 32 WCA 
particles undergoing shear flow. The crosses correspond to the directly determined probability ratio, and 



xy^ 



xy^ 



the open circles correspond to ( e 



PP VpF^-t 



pp!:(t)=A 



The inset shows the relationship between PP^^ (x) 



and PP (t) for a randomly selected set of transient trajectories, giving an indication of the relative 
direct and cross contributions to the probability ratio. The line shows the results expected if 
pP^y(T) = PPxyCt) (perfect correlation). 



the shear stress is odd under time reversal, the entropy production observed on a 
trajectory segment is only weakly correlated with the sign of the time averaged kinetic 
shear stress. In fact, quite frequently, thermodynamically favorable trajectories o(x) > 
have a kinetic contribution to the shear stress that has the opposite sign to that required to 
make o(x) > 0. The largest contribution to Pr((3P^ (x)) > actually occurs from region 



xy ' 



III. 



In Figure 3 we examine GFT for the "uncorrelated" phase function (|)(x) = V p^i(x) 



i=i 



As predicted, a phase function that is uncorrelated with the entropy production will have 



Pr((^(x) = A)/Pr((^(x) = -A) = (e"^^'^)_^^^^^=l VA. Clearly for the values of A 

tested, GFT holds, and an examination of randomly selected trajectories as shown in the 
inset confirms that most trajectories reside in regions II and III. 

We also consider two functions whose behaviour is expected to be quite different to 
that of o(x) : (^(x) = Pxy(x) = - j^dsP^yCs) and the phase function (^(x) = (pP^j,(x)) = 

1 r V 

— ds(3P^y(s) . These functions may be considered as nonlinear functions of the 
instantaneous entropy production. 
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FIGURE 3. A plot of Pr(^p^(T) = A)/Pr(^p',(T) = -A) versus A for a system of 32 WCA 

i=l i=l 

particles undergoing shear flow. The crosses correspond to the directly determined probability ratio, and 



the open circles correspond to ( e 



PPxyVFct 



E(pii)(^)=A 



The inset shows the relationship between PP^ (x) 



and V' Pxi('^) fo'" ^ randomly selected set of transient trajectories, giving an indication of the relative 



direct and cross contributions to the probability ratio. 



Equation (3) is tested for these functions in Figures 4 and 5, respectively and the relation 
is verified by the numerical data. In contrast to the previous FT (1-10), the probability 
ratio is not an exponential function of the phase function considered. 

For (|)(x) = (3P (x) it is possible to simplify Equation (3) further, obtaining. 
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FIGURE 4. A plot of Pr((PP;^^(T))= A)/Pr((PP;^^(T)) = -A) versus A for a system of 32 WCA 
particles undergoing shear flow. The crosses correspond to the directly determined probability ratio, and 
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FIGURE 5. A plot of Pr((PP^y(T))' = A)/Pr((PP^y(T))' = -A) versus A for a system of 32 WCA 
particles undergoing shear flow. The crosses correspond to the directly determined probability ratio, and 



the open circles correspond to ( e 
probability ratio by Equation (13). 



PPxyVFcT 



(PPxv(t)) =A 



The line shows the behaviour predicted for the 



P((pPxy(X)f =A) 

p((|3P,y(x)) =-A) iP'^>^'^i 



= expCpP, rx)VYT) 



1/3 \"i 



= ^exp(A^'"VYT)^ (13) 

= expC-A^^^Vyr) 

The curve shown on Figure 5 is the behaviour expected if Equation (13) is true, and the 
figure indicates that the numerical data is well determined by this relation. This function 

is not exponential in A, but is exponential in A"\ and this behaviour will persist for all 

finite x. 

We have obtained a Generalised Fluctuation Formula which is applicable to any 
phase function that is odd under a time-reversal mapping. The numerical results verify 
our prediction. Since the formula is quite general it may readily be extended to such 
problems as those considering a phase variable determined by a restricted set of 
particles, which results to a local formula (12) or to the average value of a variable over 
part of a trajectory segment, an approach that may be taken to obtain a steady state 
fluctuation formula (7,13). Finally, we note that these formulae can also be applied 
asymptotically to steady state trajectory segments rather than the transient segments 
considered here. 
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